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We investigate Lee- Yang zeros of generating functions of dynamical observables and establish a 
general relation between phase transitions in ensembles of trajectories of stochastic many-body sys- 
tems and the time evolution of high-order cumulants of such observables. This connects dynamical 
free-energies for full counting statistics in the long-time limit, which can be obtained via large- 
deviation methods and whose singularities indicate dynamical phase transitions, to observables that 
are directly accessible in simulation and experiment. As an illustration we consider facilitated spin 
models of glasses and show that from the short-time behavior of high-order cumulants it is possible 
to infer the existence and location of dynamical or "space-time" transitions in these systems. 

PACS numbers: 05.40.a, 64.70.Pf, 72.70.+m 



Introduction. — Phase transitions are a central topic 
in the statistical mechanics of equilibrium and non- 
equilibrium systems. In problems with physically mean- 
ingful interactions, phase transitions occur in the limit 
of large system size. For dynamical phase transitions, 
this also implies the limit of long times. In experiment 
or simulation of systems with complex dynamics, how- 
ever, often only the short-time dynamics can be probed, 
making it difficult to investigate dynamical transitions. 
Furthermore, such non-equilibrium transitions may be 
driven by "counting" fields [THE] which can be hard to 
relate to physically accessible parameters. In this Let- 
ter we provide a potential resolution to these problems 
by establishing a connection between phase transitions in 
ensembles of long-time dynamical trajectories of classical 
stochastic many-body systems [IMS] and the dynamics of 
physical observables at short times [9TfTT]. 

Figure 1 illustrates our approach and results. Panel (a) 
shows a dynamical trajectory of a simple lattice sys- 
tem which displays complex dynamics, in this exam- 
ple the one-dimensional East model of a glass former 
|12j . Facilitated models such as the East model show 
pronounced dynamical spatial fluctuations [13] (a phe- 
nomenon characteristic of glasses known as dynamical 
heterogeneity; for reviews see Refs. [T4HT6] ). These large 
spatio-temporal fluctuations give rise to fat tails [17] in 
the full counting statistics (FCS) [H] of time-extensive 
dynamical observables. This is shown in Fig. 1(b) for 
the dynamical activity k = K/t per unit time of the East 
model. The dynamical activity K is the number of con- 
figuration changes in a trajectory [U E] 01 [19]. Associ- 
ated with the distribution P(K, t) is the moment gener- 
ating function (MGF) Z{s,t) = e - sK P(K, t), which 
at long times t — > oo has a large-deviation (LD) form, 
Z{s,t) oc exp{i6»(s)} [IHS]. The LD function -6{s) is a 
dynamical free-energy for the counting process. Its an- 
alytic properties carry information about the phase be- 
havior of ensembles of trajectories. 
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FIG. 1. (color online) . Dynamical phase transitions and Lee- 
Yang zeros, (a) Trajectory of the one-dimensional East model, 
showing the state of up/down (black/white) spins on the lat- 
tice in time (with temperature T = 0.8, N = 150 lattices 
sites, and i max = 1000 time steps); dynamic heterogeneity 
is evident in the "space-time bubbles" of the trajectory [13] . 
(b) Probability P(k, t) of the activity k = K/t per unit time 
(green/full curve); dashed line is a Gaussian distribution with 
same mean and variance, (c) The LD function 9(s) is singular 
at s c = (blue/top) where the average activity is discontinu- 
ous (red/bottom), indicative of a first-order dynamical tran- 
sition [2]. (d) Lee- Yang zeros of the MGF in the complex-s 
plane, extracted from the cumulants of K, allow us to extrap- 
olate s c from short-time observables. 



In the East model example 6(s) has a first-order singu- 
larity at s c = 0, Fig. 1(c), which indicates that dynamics 
takes place at the coexistence of two dynamical or "space- 
time" phases, an active phase with (k)(s) = —d s 9(s) > 
for t — > oo (the equilibrium phase where relaxation is 
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possible) and an inactive phase with (k) (s) = (the dy- 
namical "glass" phase) [2j. The variable s driving the 
transition is a "counting" field which biases the trajec- 
tory ensemble from the actual dynamical one at s = 0, 
but whose connection to physically controllable param- 
eters can be hard to establish. Similar trajectory phase 
transitions are observed in other classical and quantum 
systems with complex dynamics [3HH] . 

Here we demonstrate that it is possible to infer the ex- 
istence and location of singularities of 0(s), indicative of 
phase transitions in the space of long-time trajectories, 
from short-time observables at s = 0. Specifically, we 
show that: (i) from a dynamical version of the Lee- Yang 
theorem [20] , zeros of the MGF in the complex- s plane at 
finite t will move to the real-s line in the limit of t — > oo 
if there are any singularities in #(s); and (ii) these zeros 
can be obtained from the short-time and finite-size be- 
havior of cumulants [9HTT] of dynamic observables such 
as the activity. Figure 1(d) illustrates this result for the 
East model: the s c = singularity of the thermodynamic 
and long-time limit can be extrapolated from the leading 
Lee- Yang zeros extracted from short-time cumulant dy- 
namics. This offers the possibility of studying trajectory 
phase transitions in FCS via observables that are directly 
accessible in simulation and experiment. 

Formalism. — For concreteness we consider stochastic 
processes described by the Master equation [2T] 

d t P(C, t) = -r(C)P(C, t) + J2 w ( c ' -> Q p (C, t). (1) 

c 

Here, P(C,t) is the probability that the system is in the 
configuration C at time t. The transition rate from con- 
figuration C to C is denoted as W(C -> C) and r(C) = 
J2c> W(C — > C) is the total escape rate from C. By def- 
inition W(C — > C) = 0. Equation ([I]) can be written 
in the convenient matrix notation dt\P(t)) — W\P(t)), 
where the matrix W is defined as 

W(C,C') = W(C ^C)-r(C)8c,c, (2) 

and the vector \P(t)} contains the probabilities P(C,i)'s. 

We classify trajectories according to their dynamical 
activity K — the total number of spin-flips in the case 
of spin models considered here [H [2J. (Similar argu- 
ments can be applied to analyze ensembles of trajecto- 
ries classified by other time-extensive dynamic observ- 
ables, see e. g. Refs. [IH3])- The probability that the 
system is in configuration C at time t, having changed 
configuration K times, is denoted as P(C\K,t). Then 
P(K,t) = Ec P (C\K,t) and Z(s,t) = T,c P {^s,t), 
where P(C,s,t) = J2k p ( C \ k , t)e~ sK [HE]. The corre- 
sponding vector \P(s,t)) obeys d t \P(s,t)) =W s \P(s,t)), 
where the generalized Master operator is [TJ [2] 

W.(C,C) = e- s W(C -> C) - r(C)S c ,a- (3) 

Formally, the solution to Eq. Q is \P(s, f)) = e w =*|P(0)), 
assuming for instance that the initial state \P(0)} is the 



equilibrium distribution defined by W s= o|-P(0)) = 0. By 
using the "flat" state, (— | = (1,...,1), we can express 
the MGF as Z(s,t) = (-\P(s,t)) = (-|e w '*|P(0)) = 
J2j -,Cj(s)e Xi ^ t in terms of the eigenvalues Xj{s) of W s 
and corresponding expansion coefficients Cj(s). The cu- 
mulant generating function (CGF) is defined in terms 
of the MGF as Q(s,t) = log Z(s,t), which delivers the 
cumulants of K by differentiation with respect to the 
counting variable s at s = 0, 

«x»»(*) = (-i)»a s "e( S) i)|^ . (4) 

At long times the MGF function becomes exponential in 
time pQ ; its rate of change is determined by the eigenvalue 
with the largest real-part, such that Q(s,t) — > t9(s), 
where 6(s) = max[A :) (s)] is the LD function. 

Singularities and dynamical transitions. — Fluctua- 
tions in the dynamical system can be understood from 
the analytic properties of 9(s). For example, a first- 
order dynamical phase transition corresponds to singu- 
larities in 9(s) so that its first derivative is discontinuous 
[2], see Fig. 1(c). This occurs at a real s — s c where 
the two largest eigenvalues of W s become degenerate, 
^o( s c) = ^i(sc)- As a central result of this work, we 
show below how such dynamical phase transitions, oc- 
curring in the long-time limit, can be inferred from the 
high-order cumulants of K at finite times and at s = 0, 
i. e. evolving under the unbiased dynamics. 

To this end we consider the zeros of the MGF in the 
vicinity of the transition value, s ~ s Cl where the two 
largest eigenvalues are nearly degenerate Xq(s) ~ Ai(s) 
and we may write Z(s,t) ~ c (s)e Ao(s)t + c 1 (s)e Xl ^ t . 
The zeros of the MCF are determined by the equations 
Ao(s) = Ai(s) + [logci(s)/c (s)-(-i7r(2TO-|-l)]/i for integer 
m. In the long-time limit, these equations all reduce to 
Ao(s) = Ai(s), and thus with increasing time all zeros 
Sj(t) move towards the transition value s c on the real- 
axis. (At finite times, the zeros must be complex, since 
Z(s,t) > for real s.) This is in essence the theory 
of phase transitions of Lee and Yang [2D], here applied 
to dynamical systems [22]. Accordingly, we refer to the 
(time-dependent) zeros Sj(t) of the MGF as Lee- Yang 
zeros. 

High- order cumulants and Lee-Yang zeros. — The mo- 
tion of the Lee- Yang zeros in the complex plane can be 
inferred from the high-order cumulants of K. Impor- 
tantly, the zeros of the MGF correspond to logarithmic 
singularities of the CGF which determine the high-order 
derivatives of the CGF (the cumulants) according to Dar- 
boux's theorem [23]. Writing the MGF in terms of the 
Lee- Yang zeros as Z(s,t) = Ylj[sj(t) — s]/sj(t), where 
Z(0, t) = 1 reflects the normalization J^k p (Kj t) = 1 
at all times, the CGF becomes 6(s,i) = J2j{\og{sj(t) — 
s] — \og[sj(t)]). The Lee-Yang zeros come in complex- 
conjugate pairs, since the MGF is real for real s. Com- 



bined with Eq. Q we readily find [SHIT] 

-marg[sj(t)] 

««*»(*) = (-l^n-lllg |gj(t)|n ■ (5) 

This result shows that higher-order cumulants generically 
grow as the factorial of the cumulant order n, and os- 
cillate as a function of any parameter that changes the 
complex argument arg{ Sj(i)} This behavior has been 
observed experimentally [SJ[53]. For large n, the sum is 
dominated by the pair So(t) and Sg(t) of zeros closest to 
s = 0, and the expression further simplifies to [9HTT] [25] 
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(CT(^(-l)(^)(n-l)! ^^M . (6) 

We can solve this simple relation for sq, given the ratios 
of cumulants dt\t) = ((K n±1 )){t)/((K n ))(t). We then 
obtain the matrix equation 
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which directly yields so(t) from four consecutive cumu- 
lants [TO] [TTJ [26]. We now employ this method to in- 
vestigate dynamical phase transitions in kinetically con- 
strained models of glass formers. 

Dynamical transitions in facilitated glass models. — 
As an example of how the ideas above can be applied, 
we study trajectory transitions [S] in facilitated spin 
models of glasses [12]. For simplicity we consider one- 
dimensional models, defined in terms of binary variables 
rij = 0, 1, where i = 1, . . . , N denote sites on a chain. 
The energy function is E = JJ2i n ii an< ^ a U interactions 
emerge via kinetic constrains, which stipulate that a site 
i changes with a rate that is determined by the state of 
its nearest neighbors i ± 1 [12] . Concretely, we focus on 
the Fredrickson-Andersen (FA) model [37] and on the 
East model [35J . In the FA model a site can only change 
if either of its nearest neighbors is in the up state, i. e. 
the transitions 11 — > 10 and 11 — > 01 occur with rate 1, 
11 <- 10 and 11 <- 01 with rate e~ J / T ', but 010 <=> 000 
are not allowed. In the East model facilitation is via the 
left neighbor only, so that 11 — s- 10 and 11 <— 10 occur 
with rates 1 and e~ J l T , respectively, but 010 ^ 000 and 
011 ^ 001 are not allowed. At low T, there is a conflict 
between lowering the energy and having enough excited 
spins to evolve dynamically, which gives rise to glassy 
slow-down and dynamical heterogeneity [T5J [T3] in these 
systems; the East model in particular seems to capture 
the basic physics of glassy dynamical arrest [TS] . 

Results. — Figure [2] shows our numerical simulations 
for the high-order cumulants of the activity K as func- 
tions of time for the East and FA models (full lines). 
The cumulants grow dramatically with the cumulant or- 
der and oscillate as functions of time (the absolute value 
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FIG. 2. (color online). High-order cumulants of the activity in 
facilitated models of glass formers. The upper (lower) panels 
show the time-evolution of the cumulants of order n = 4 — 7 
for the East (FA) model with TV = 30 sites at T = 0.8 (left) 
and T = 1.0 (right). Simulations are shown with full lines, 
while dashed lines indicate the approximation Eq. |6| based 
on the closest pair of Lee- Yang zeros which are extracted from 
the numerical data using Eq. |7b. The motion of the closest 
pair of Lee- Yang zeros corresponding to the East model at 
T = 0.8 is shown in Fig. [jjd). 



is shown on a logarithmic scale, such that downwards- 
pointing spikes on the curves correspond to the cumu- 
lants crossing zero). This is due to the Lee- Yang zeros 
Sj(t) approaching the transition value at s c = accord- 
ing to Eq. g), causing the large growth of the cumulants. 
Initially, P(K, t — 0) = Sk.o an d all cumulants of the 
activity are zero, implying that the Lee- Yang zeros are 
infinitely far from s c = and l/\sj(t = 0)| = 0. At very 
short time, where P{K = 0, t) « 1 > P(l, t) > P(2, t) > 
. . ., the leading pair of Lee- Yang zeros are determined by 
the equation Z(s, t) ~ P(0, t) + P(l, t)e~ s = with solu- 
tions s (f), s (f) =-log{P(0,t)/P(l,i)}±i7r. Thus, to 
begin with the Lee- Yang zeros move along the lines ±itt 
from — oo ± i7r, before approaching s c = 0. We now use 
Eq. fff\ to deduce the motion of the leading Lee- Yang 
zeros from the numerical data. 

Figure [ljd) shows the leading pair of Lee- Yang ze- 
ros, so(t) and Sg(t), for the East model with N = 30 
sites at temperature T — 0.8 as they move towards the 
first-order transition point at s c = 0. To validate the 
extraction of the leading Lee- Yang zeros from the cumu- 
lants of the activity using Eq. Q , we plug the solution 
So(t) back into Eq. ^ and compare the result with the 
numerical data. In Fig. [2j we show the numerical results 
(full lines) together with the approximation in Eq. ^ 
based on the extracted pair of Lee- Yang zeros (dashed 
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FIG. 3. (color online). Finite-size scaling of the first-order 
transition point. Upper (lower) panels show the real-part of 
s c (t), extracted from the cumulants for the East (FA) model 
with N = 15 - 50 sites at T = 0.8 (left) and T = 1.0 (right). 
Numerical values are shown with circles, while linear fits based 
on Eq. (JsJ) are shown with full lines. Upper insets show the 
imaginary part of s c (t) as a function of the inverse time 1/i. 
Lower insets show the extrapolated long-time limit si°°\N) 
as a function of the inverse system size 1/N. 



extrapolated values of (N) as a function of the in- 
verse system size 1/N. These results show that the value 
s c = is approached in the large-system-size limit. Some 
deviations are seen for the larger systems as we reach the 
limits of the numerical accuracy of our method. Our re- 
sults show that it is possible to infer the existence and 
location of dynamical singular points, which are indica- 
tive of phase transitions in the space of long-time trajec- 
tories, from high-order short-time cumulants at s = 0. 
Our method can also be extended to systems where the 
transition point on the real-s line is at s c ^ |S]. 

Conclusions. — We have investigated the Lee- Yang ze- 
ros of generating functions of dynamical observables and 
demonstrated how singularities in the long-time limit, in- 
dicative of dynamical phase transitions, can be inferred 
from the short-time dynamics of high-order cumulants 
in finite-size systems. We hope that our approach may 
facilitate theoretical and experimental studies of trajec- 
tory phase transitions in stochastic many-body systems. 
An important task to address in future work is to apply 
similar ideas to dynamical phase transitions in quantum 
many-body systems pj [5] • 
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line). The figure corroborates that we indeed are extract- 
ing the leading pair of Lee- Yang zeros. Some deviations, 
in particular at long times, are observed as the second 
pair of Lee- Yang zeros also come close to s = and start 
contributing significantly to the sum in Eq. jBJ. bince 
the second pair of Lee- Yang zeros is not included in Eq. 

a shift in the frequency of the oscillations as a func- 
tion of time is also observed. If needed, the accuracy of 
the method can be improved by using higher cumulants 

union lag . 

In Fig. [3] we analyze the finite-size scaling of the tran- 
sition point s c . The real-part of the transition point is 
predicted to scale as 



Re[s c (t)] -a/t + s^iN), 



(8) 



where the coefficient a depends on the temperature T, 
and si°°\N) oc 1/N is the long-time value, which for 
the East and FA models should approach s c = in the 
limit N — > 00. Our numerical results for the East and 
FA models confirm the predicted scaling behavior. For 
each system size in the range N = 15 to 50, we find 
an approximately linear dependence on the inverse time 
1/t, allowing us to extrapolate the values of si^^N) in 
the t — ¥ 00 limit. We also verify that the imaginary 
part of s c approaches zero in the long-time/large-system 
limit, see upper insets. In the lower insets, we show the 
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